home *** CD-ROM | disk | FTP | other *** search
/ Turnbull China Bikeride / Turnbull China Bikeride - Disc 1.iso / ARGONET / PD / MATHS / RLAB / RLAB125.ZIP / !RLaB / rlib / plplot < prev    next >
Text File  |  1996-03-01  |  58KB  |  2,602 lines

  1. #
  2. # New plot.r for use with PLPLOT library.
  3. # The help files for these functions are in
  4. # misc/plhelp
  5. #
  6.  
  7. # plplot.r
  8.  
  9. # This file is a part of RLaB ("Our"-LaB)
  10. # Copyright (C) 1994  Ian R. Searle
  11.  
  12. # This program is free software; you can redistribute it and/or modify
  13. # it under the terms of the GNU General Public License as published by
  14. # the Free Software Foundation; either version 2 of the License, or
  15. # (at your option) any later version.
  16.  
  17. # This program is distributed in the hope that it will be useful,
  18. # but WITHOUT ANY WARRANTY; without even the implied warranty of
  19. # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  20. # GNU General Public License for more details.
  21.  
  22. # You should have received a copy of the GNU General Public License
  23. # along with this program; if not, write to the Free Software
  24. # Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
  25.  
  26. # See the file ../COPYING
  27.  
  28. #
  29. # If your system does not deal with Infs and NaNs well, then
  30. # uncomment the following lines.
  31. #
  32.  
  33. # static (isinf, isnan)
  34. #
  35. # isinf = function ( A ) { return (0); };
  36. # isnan = function ( A ) { return (0); };
  37.  
  38. static (WIN)        # The static plot window structure
  39. static (P)        # The active/current plot window
  40.  
  41. #
  42. # Maintain the transformations for 3-D plots.
  43. #
  44.  
  45. static (basex, basey, height)
  46. basex = 2;
  47. basey = 2;
  48. height = 4;
  49.  
  50. #
  51. # Static (private) functions. For use from within
  52. # this file only.
  53. #
  54.  
  55. static (create_plot_object)
  56. static (check_plot_object)
  57. static (xy_scales)
  58. static (x_scales)
  59. static (y_scales)
  60. static (z_scales)
  61. static (XYZ_scales)
  62. static (list_scales)
  63. static (list_sort)
  64. static (hist_scales)
  65. static (plot_matrix)
  66. static (plot_list)
  67. static (check_3d_list)
  68. static (find_char)
  69. static (get_style)
  70. static (make_legend)
  71. static (plhold_first)
  72.  
  73. #
  74. # Defaults
  75. #
  76.  
  77. static (grid_x_default, grid_y_default)
  78. static (grid_3x_default, grid_3y_default, grid_3z_default)
  79.  
  80. grid_x_default = "bcnst";
  81. grid_y_default = "bcnstv";
  82. grid_3x_default = "bnstu";
  83. grid_3y_default = "bnstu";
  84. grid_3z_default = "bcdmnstuv";
  85.  
  86. static (subplot_f)
  87. subplot_f = 0;
  88.  
  89. #
  90. # Create the default plot-object.
  91. # Initialize to all the default values
  92. #
  93.  
  94. if (!exist (WIN)) 
  95. {
  96.   # Create the plot-object list
  97.   WIN = <<>>;
  98. }
  99.  
  100. create_plot_object = function ( N, nx, ny )
  101. {
  102.   if (!exist (N)) { N = 0; }
  103.   
  104.   pobj = <<>>;
  105.   pobj.subplot = 0;        # The current subplot no.
  106.   pobj.nplot = nx*ny;        # Total no. of plots on window
  107.   
  108.   pobj.fontld = 0;        # Loaded extended fonts?
  109.   
  110.   for (i in 1:(nx*ny))
  111.   {
  112.     pobj.style.[i] = "line";    # The type/style of plot to draw
  113.     pobj.nbin[i] = 1j;          # The number of bins for a histogram
  114.     pobj.width[i] = 1;        # The pen width for current plot
  115.     pobj.font[i] = 1;        # The current font
  116.     pobj.xlabel[i] = "";
  117.     pobj.ylabel[i] = "";
  118.     pobj.zlabel[i] = "";
  119.     pobj.title[i] = "";
  120.     pobj.orientation[i] = "portrait";
  121.     pobj.desc.[i] = "default";        # The legend description
  122.     pobj.gridx[i] =  grid_x_default;    # Plot axes style, 2D-X
  123.     pobj.gridy[i] =  grid_y_default;    # Plot axes style, 2D-Y
  124.     pobj.grid3x[i] = grid_3x_default;    # Plot axes style, 3D-X
  125.     pobj.grid3y[i] = grid_3y_default;    # Plot axes style, 3D-Y
  126.     pobj.grid3z[i] = grid_3z_default;    # Plot axes style, 3D-Z
  127.     pobj.aspect[i] = 0;                # Plot aspect style
  128.     pobj.alt[i] = 60;
  129.     pobj.az[i] = 45;
  130.     
  131.     pobj.xmin[i] = 1j;
  132.     pobj.xmax[i] = 1j;
  133.     pobj.ymin[i] = 1j;
  134.     pobj.ymax[i] = 1j;
  135.     pobj.zmin[i] = 1j;
  136.     pobj.zmax[i] = 1j;
  137.     
  138.     pobj.page.xp = 0;
  139.     pobj.page.yp = 0;
  140.     pobj.xleng = 400;
  141.     pobj.yleng = 300;
  142.     pobj.xoff = 200;
  143.     pobj.yoff = 200;
  144.  
  145.     pobj.color[i;] = 1:14;               # 14 possible colors...
  146.     pobj.lstyle[i;] = 1:8;               # 8 possible line styles...
  147.     pobj.pstyle[i;] = 1:8;               # 8 possible point styles...
  148.   }
  149.   
  150.   #
  151.   # Save the newly generated plot-object
  152.   # in a list of plot-objects.
  153.   
  154.   WIN.[N] = pobj;
  155. };
  156.  
  157. ##############################################################################
  158. #
  159. # Check to make sure a plot-object exists. If one
  160. # does not exist, create it.
  161. #
  162.  
  163. check_plot_object = function ()
  164. {
  165.   if (length (WIN) == 0)
  166.   {
  167.     plstart();
  168.     return 0;
  169.   }
  170.   return 1;
  171. };
  172.  
  173. ##############################################################################
  174. #
  175. # Set the current plot window
  176. # Default value = 0
  177. #
  178.  
  179. plwin = function ( N )
  180. {
  181.   check_plot_object ();
  182.   if (!exist (N)) { N = 0; }
  183.   
  184.   # Check to make sure N is valid
  185.   
  186.   for (i in members (WIN))
  187.   {
  188.     if (N == strtod (i))
  189.     {
  190.       _plsstrm (N);
  191.       return P = N;
  192.     }
  193.   }
  194.   printf ("plwin: invalid argument, N = %i\n", N);
  195.   printf ("      valid values are:\n");
  196.   WIN?
  197. };
  198.  
  199. ##############################################################################
  200. #
  201. # Show the current plot-window, and the possibilities
  202. #
  203.  
  204. showplwin = function ( all )
  205. {
  206.   if (length (WIN) == 0)
  207.   {
  208.     printf ("No plot objects\n");
  209.     return 0;
  210.   }
  211.   
  212.   printf ("Current plot-window is:\t\t%i\n", P);
  213.   printf ("Available plot windows are:\t");
  214.   for (i in members (WIN))
  215.   {
  216.     printf ("%s   ", i);
  217.   }
  218.   printf ("\n");
  219.  
  220.   if (exist (all))
  221.   {
  222.     for (i in members (WIN.[P]))
  223.     {
  224.       WIN.[P].[i]
  225.     }
  226.   }
  227. };
  228.  
  229. getplot = function ( win_no )
  230. {
  231.   local (win_no)
  232.  
  233.   if (length (WIN) != 0)
  234.   {
  235.     if (!exist (win_no)) { win_no = P; }
  236.  
  237.     if (exist (WIN.[win_no]))
  238.     {
  239.       return (WIN.[win_no]);
  240.     else
  241.       return 0;
  242.     }
  243.   }
  244.   return <<>>;
  245. };
  246.  
  247.  
  248. ##############################################################################
  249. #
  250. # Set/start/select the plot device
  251. #
  252.  
  253. plstart = function ( nx, ny, dev )
  254. {
  255.   if (!exist (nx)) { nx = 1; }
  256.   if (!exist (ny)) { ny = 1; }
  257.   if (!exist (dev)) { dev = "?"; }
  258.   
  259.   # Create the plot-object
  260.   # First, figure out the index
  261.   if (!exist (P))
  262.   {
  263.     P = 0;
  264.   else
  265.     P = P + 1;
  266.   }
  267.   
  268.   create_plot_object (P, nx, ny);
  269.   _plsstrm (P);
  270.   
  271.   _plscolbg (0,0,0); # white
  272.   # Default window size for X
  273.   _plspage (0, 0, 400, 300, 200, 200);
  274.   
  275.   # Start up the plot-window
  276.   _plstart (dev, nx, ny);
  277.   
  278.   _plwid (8);
  279.   
  280.   # Turn between plot pause off
  281.   _plspause (0);
  282.   _pltext ();
  283.   
  284.   return P;
  285. };
  286.  
  287. ##############################################################################
  288. #
  289. # Close a plot device. We must destroy the current plot-object
  290. # And switch the output stream back to the default.
  291. #
  292.  
  293. plclose = function ()
  294. {
  295.   if (size (WIN) > 1)
  296.   {   
  297.     #
  298.     # Clear WIN.[P] and reset P to 1st plot-window
  299.     #
  300.  
  301.     clear (WIN.[P]);
  302.     _plend1 ();
  303.     _plsstrm (strtod (members (WIN)[1]));
  304.     P = strtod (members (WIN)[1]);
  305.     return P;
  306.  
  307.   else if (size (WIN) == 1) {
  308.  
  309.     if (exist (WIN.[P])) 
  310.     { 
  311.       clear (WIN.[P]); 
  312.       clear (P);
  313.     }
  314.     _plend1 ();
  315.     return 1;
  316.  
  317.   else if (size (WIN) == 0) {
  318.  
  319.     return 0;
  320.  
  321.   } } }
  322. };
  323.  
  324. ##############################################################################
  325. #
  326. # Close ALL the plot-windows
  327. #
  328.  
  329. plend = function ()
  330. {
  331.   _plend ();
  332.   if (exist (WIN)) { clear (WIN); }
  333.   if (exist (P)) { clear (P); }
  334.   WIN = <<>>;
  335. };
  336.  
  337. ##############################################################################
  338. #
  339. # Change plot aspect ratio
  340. #
  341.  
  342. plaspect = function ( aspect )
  343. {
  344.   check_plot_object ();
  345.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  346.   if (!exist (aspect))
  347.   {
  348.     WIN.[P].aspect[i] = 0;
  349.   else
  350.     WIN.[P].aspect[i] = aspect;
  351.   }
  352. };
  353.  
  354. ##############################################################################
  355. #
  356. # Change plot line style
  357. #
  358.  
  359. plstyle = function ( style )
  360. {
  361.   check_plot_object ();
  362.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  363.  
  364.   if (exist (style)) 
  365.   {
  366.     if (class (style) == "string") 
  367.     {
  368.       for (j in 1:style.n)
  369.       {
  370.     if (style[j] != "line" && ...
  371.             style[j] != "point" && ...
  372.             style[j] != "line-point")
  373.     {
  374.       error ("plstyle: STYLE must be either " + ...
  375.              "\"point\", \"line\" or \"line-point\"");
  376.     }
  377.       }
  378.       WIN.[P].style.[i] = style;
  379.     }
  380.     return 1;
  381.   }
  382.   WIN.[P].style.[i] = "line";
  383.   return 1;
  384. };
  385.  
  386. ##############################################################################
  387. #
  388. # Control of the plot line style.
  389. # There are 8 line styles
  390.  
  391. plline = function ( line_style )
  392. {
  393.   check_plot_object ();
  394.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  395.  
  396.   if (exist (line_style)) 
  397.   {
  398.     if (class (line_style) == "num") 
  399.     {
  400.       if (line_style.nc != 8)
  401.       {
  402.     error ("plpoint: LVEC must be 1x8 in size");
  403.       }
  404.       WIN.[P].lstyle[i;] = line_style;
  405.     }
  406.     return 1;
  407.   }
  408.   WIN.[P].lstyle[i;] = 1:8;
  409.   return 1;
  410. };
  411.  
  412. ##############################################################################
  413. #
  414. # Control of the plot point style.
  415. # There are 8 line styles
  416.  
  417. plpoint = function ( point_style )
  418. {
  419.   check_plot_object ();
  420.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  421.  
  422.   if (exist (point_style)) 
  423.   {
  424.     if (class (point_style) == "num") 
  425.     {
  426.       if (point_style.nc != 8)
  427.       {
  428.     error ("plpoint: PVEC must be 1x8 in size");
  429.       }
  430.       WIN.[P].pstyle[i;] = point_style;
  431.     }
  432.     return 1;
  433.   }
  434.   WIN.[P].pstyle[i;] = 1:8;
  435.   return 1;
  436. };
  437.  
  438. ##############################################################################
  439. #
  440. # Get the right value of line-style
  441. #
  442.  
  443. get_style = function ( STY, K )
  444. {
  445.   local (sty);
  446.   sty = mod(K, STY.n);
  447.   if(sty == 0) 
  448.   { 
  449.     sty = STY.n; 
  450.   }
  451.   return STY[sty];
  452. };
  453.  
  454. ##############################################################################
  455. #
  456. # Change fonts
  457. #
  458.  
  459. plfont = function ( font )
  460. {
  461.   check_plot_object ();
  462.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  463.   
  464.   if (!exist (font)) { font = 1; }
  465.   
  466.   if (WIN.[P].fontld == 0)
  467.   {
  468.     _plfontld (1);
  469.     WIN.[P].fontld = 1;
  470.   }
  471.   
  472.   WIN.[P].font[i] = font;
  473.   return P;
  474. };
  475.  
  476. ##############################################################################
  477. #
  478. # Change pen width
  479. #
  480.  
  481. plwid = function ( width )
  482. {
  483.   check_plot_object ();
  484.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  485.   
  486.   if (!exist (width)) { width = 1; }
  487.   WIN.[P].width[i] = width;
  488.   return P;
  489. };
  490.  
  491. ##############################################################################
  492. #
  493. # Place some text on the plot
  494. #
  495.  
  496. plptex = function ( text, x , y , dx , dy , just )
  497. {
  498.   if (!check_plot_object ()) 
  499.   {
  500.     printf ("Must use plot() before plptex()\n");
  501.     return 0;
  502.   }
  503.   
  504.   if (!exist (x)) { x = 0; }
  505.   if (!exist (y)) { y = 0; }
  506.   if (!exist (dx)) { dx = abs(x)+1; }
  507.   if (!exist (dy)) { dy = 0; }
  508.   if (!exist (just)) { just = 0; }
  509.   
  510.   _plptex (x, y, dx, dy, just, text);
  511. };
  512.  
  513. ##############################################################################
  514. #
  515. # Set up the viewing altitude for 3-D plots
  516. #
  517.  
  518. plalt = function ( ALT )
  519. {
  520.   check_plot_object ();
  521.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  522.  
  523.   if (exist (ALT)) 
  524.   { 
  525.     WIN.[P].alt[i] = ALT; 
  526.   else
  527.     WIN.[P].alt[i] = 60;
  528.   }
  529.   return P;
  530. };
  531.  
  532. ##############################################################################
  533. #
  534. # Set the viewing azimuth for 3-D plots
  535. #
  536.  
  537. plaz = function ( AZ )
  538. {
  539.   check_plot_object ();
  540.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  541.  
  542.   if (exist (AZ)) 
  543.   { 
  544.     WIN.[P].az[i] = AZ; 
  545.   else
  546.     WIN.[P].az[i] = 45;
  547.   }
  548.   return P;
  549. };
  550.  
  551. ##############################################################################
  552. #
  553. # Find a character in a string
  554. #
  555.  
  556. find_char = function ( str , char )
  557. {  
  558.   tmp = strsplt (str);
  559.   for (i in 1:tmp.n)
  560.   {
  561.     if (tmp[i] == char) 
  562.     {
  563.       return i;
  564.     }
  565.   }
  566.   return 0;
  567. };
  568.  
  569. ##############################################################################
  570. #
  571. # Sort list element names/labels by numeric order, then string order.
  572. #
  573.  
  574. list_sort = function ( L )
  575. {
  576.   tl = <<>>;
  577.   j = k = 1;
  578.  
  579.   for (i in members (L))
  580.   {
  581.     if (!isnan (strtod (i)))
  582.     {
  583.       num[j] = i;
  584.       j++;
  585.       else
  586.       char[k] = i;
  587.       k++;
  588.     }
  589.   }
  590.   
  591.   # Sort the numeric labels
  592.   
  593.   if (exist (num))
  594.   {
  595.     num = sort (strtod (num)).val;
  596.     tl.num = num;
  597.   }
  598.  
  599.   if (exist (char))
  600.   {
  601.     tl.char = char;
  602.   }
  603.  
  604.   return tl;
  605. };
  606.  
  607. ##############################################################################
  608. #
  609. # Set the subplot, this overides the action in plot().
  610. #
  611.  
  612. subplot = function ( sub )
  613. {
  614.   check_plot_object ();
  615.   
  616.   if (!exist (sub))
  617.   {
  618.     subplot_f = 0;
  619.     _pladv (0);
  620.   else
  621.     if (sub > WIN.[P].nplot)
  622.     {
  623.       error ("Current window does not have this many subplots");
  624.     }
  625.     if (sub > 0)
  626.     {
  627.       WIN.[P].subplot = sub - 1;
  628.       subplot_f = 1;
  629.       _pladv (sub);
  630.     else
  631.       if (sub == 0)
  632.       {
  633.         # Do not advance, stay at current subplot
  634.         WIN.[P].subplot = WIN.[P].subplot - 1;
  635.         subplot_f = 1;
  636.       }
  637.     }
  638.   }
  639. };
  640.  
  641. ##############################################################################
  642. #
  643. # Plot the columns of a matrix (X-Y plot).
  644. #
  645. ##############################################################################
  646.  
  647. plot = function ( data, key, textf )
  648. {
  649.   check_plot_object ();
  650.   
  651.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  652.   if (!exist (key)) { key = 1; }
  653.   if (!exist (textf)) { textf = 1; }
  654.  
  655.   #
  656.   # Draw the graph
  657.   # Step through the matrix plotting
  658.   # each column versus the 1st
  659.   #
  660.   
  661.   if (class (data) == "num")
  662.   {
  663.     #
  664.     # Set up the plot basics
  665.     #
  666.  
  667.     if (abs (key) > data.nc)
  668.     {
  669.       error_1 ("plot: KEY argument > M.nc");
  670.     }
  671.     
  672.     _plgra ();
  673.     _plcol (1);
  674.     _pllsty (1);
  675.     _plfont (WIN.[P].font[p]);
  676.     _plwid (WIN.[P].width[p]);
  677.  
  678.     if (!subplot_f) 
  679.     {
  680.       _pladv (0);        # Advance 1 subplot
  681.     else
  682.       subplot_f = 0;     # The user has set the subplot
  683.     }
  684.  
  685.     if (WIN.[P].aspect[p] != 0)
  686.     {
  687.       _plvasp (WIN.[P].aspect[p]);
  688.     else
  689.       _plvsta ();
  690.     }
  691.     
  692.     #
  693.     # Compute scale limits
  694.     #
  695.  
  696.     k = find ((1:data.nc) != abs (key));
  697.     if (key > 0)
  698.     {
  699.       if (data.nc != 1)
  700.       {
  701.     x_scales ( real(data)[;key], p, xmin, xmax );
  702.     y_scales ( real(data)[;k],   p, ymin, ymax );
  703.       else
  704.     x_scales ( (1:data.nr)', p, xmin, xmax );
  705.     y_scales ( real(data),   p, ymin, ymax );
  706.       }
  707.     else if (key < 0) {
  708.       x_scales ( real(data)[;k],   p, xmin, xmax );
  709.       y_scales ( real(data)[;abs(key)], p, ymin, ymax );
  710.     else
  711.       x_scales ( (1:data.nr)', p, xmin, xmax );
  712.       y_scales ( real(data),   p, ymin, ymax );
  713.     } }
  714.  
  715.     _plwind (xmin, xmax, ymin, ymax);
  716.     _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  717.     if (plot_matrix ( data, key, p, 0, xmin, xmax, ymin, ymax, ymax-ymin ) < 0) 
  718.     {
  719.       return -1;
  720.     }
  721.     
  722.     else if (class (data) == "list") {
  723.       
  724.       _plgra ();
  725.       _plcol (1);
  726.       _pllsty (1);
  727.       _plfont (WIN.[P].font[p]);
  728.       _plwid (WIN.[P].width[p]);
  729.  
  730.       list_scales ( data, key, p, xmin, xmax, ymin, ymax );
  731.  
  732.       if (!subplot_f) 
  733.       {
  734.     _pladv (0);        # Advance 1 subplot
  735.       else
  736.     subplot_f = 0;     # The user has set the subplot
  737.       }
  738.  
  739.       if (WIN.[P].aspect[p] != 0)
  740.       {
  741.     _plvasp (WIN.[P].aspect[p]);
  742.       else
  743.        _plvsta ();
  744.       }
  745.       
  746.       _plwind (xmin, xmax, ymin, ymax);
  747.       _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  748.       if (plot_list ( data, key, p, xmin, xmax, ymin, ymax ) < 0) 
  749.       { 
  750.     return -1;
  751.       }
  752.  
  753.     else
  754.       error ("plot: un-acceptable argument");
  755.   } }
  756.  
  757.   _pllsty (1);
  758.   _plcol (1);
  759.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  760.   _plflush ();
  761.   if (textf)
  762.   {
  763.     _pltext ();
  764.   }
  765.   
  766.   #
  767.   # Increment the plot no. so that next time
  768.   # we use the correct settings.
  769.   #
  770.   
  771.   WIN.[P].subplot = WIN.[P].subplot + 1;
  772.   return P;
  773. };
  774.  
  775. #
  776. # plhold:
  777. # Plot some data, and "hold" on for more.
  778. # Plot the data, setting up the plot as usual the first time.
  779. # On subsequent plots do not do any setup, just plot some
  780. # more. plhold_off must be called to finish up.
  781. #
  782.  
  783. plhold_first = 1;    # True (1) if plhold() has NOT been used.
  784.                         # Or if plhold_off() has been used.
  785.             # False (0) if plhold is in use
  786.  
  787. static (hxmin, hxmax, hymin, hymax)
  788.  
  789. plhold = function ( data, key )
  790. {
  791.   check_plot_object ();
  792.   
  793.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  794.   if (!exist (key)) { key = 1; }
  795.   if (abs (key) > data.nc)
  796.   {
  797.     error_1 ("plot: KEY argument > M.nc");
  798.   }
  799.   
  800.   if (plhold_first)
  801.   {
  802.     if (class (data) == "num")
  803.     {
  804.       #
  805.       # Do the setup ONCE
  806.       #
  807.       hxmin = hxmax = hymin = hymax = 0;
  808.       _plgra ();
  809.       _plcol (1);
  810.       _pllsty (1);
  811.       _plfont (WIN.[P].font[p]);
  812.       _plwid (WIN.[P].width[p]);
  813.  
  814.       if (!subplot_f) 
  815.       {
  816.     _pladv (0);        # Advance 1 subplot
  817.       else
  818.     subplot_f = 0;     # The user has set the subplot
  819.       }
  820.  
  821.       if (WIN.[P].aspect[p] != 0)
  822.       {
  823.     _plvasp (WIN.[P].aspect[p]);
  824.       else
  825.     _plvsta ();
  826.       }
  827.  
  828.       xy_scales ( real(data), p, hxmin, hxmax, hymin, hymax );
  829.       
  830.       k = find ((1:data.nc) != abs (key));
  831.       if (key > 0)
  832.       {
  833.     x_scales ( real(data)[;key], p, xmin, xmax );
  834.     if (data.nc != 1)
  835.     {
  836.       y_scales ( real(data)[;k],   p, ymin, ymax );
  837.         else
  838.       y_scales ( (1:data.nr)',   p, ymin, ymax );
  839.     }
  840.       else if (key < 0) {
  841.     x_scales ( real(data)[;k],   p, xmin, xmax );
  842.     y_scales ( real(data)[;abs(key)], p, ymin, ymax );
  843.       else
  844.     x_scales ( (1:data.nr)', p, xmin, xmax );
  845.     y_scales ( real(data),   p, ymin, ymax );
  846.       } }
  847.  
  848.       _plwind (hxmin, hxmax, hymin, hymax);
  849.       _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  850.       _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);
  851.     else
  852.       error ("plot: un-acceptable argument");
  853.     }
  854.     plhold_first = 0;
  855.   }
  856.   
  857.   if (plot_matrix ( data, key, p, 0, hxmin, hxmax, hymin, hymax, hymax-hymin ) < 0) 
  858.   { 
  859.     return -1; 
  860.   }
  861.  
  862.   _plcol (1);
  863.   _plflush ();
  864.   _pltext ();
  865.  
  866.   return P;
  867. };
  868.  
  869. ##############################################################################
  870. #
  871. # Clean up the plotting environment and get ready
  872. # for normal interactive usage.
  873. #
  874.  
  875. plhold_off = function ( )
  876. {
  877.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  878.   plhold_first = 1;  
  879.   _plcol (1);
  880.   _plflush ();
  881.   _pltext ();
  882.   WIN.[P].subplot = WIN.[P].subplot + 1;
  883.   return P;
  884. };
  885.  
  886. ##############################################################################
  887. #
  888. # Plot a 3-Dimensional graph. The data is composed in a list, with
  889. # elements `x', `y', and `z'. x and y are single-dimension arrays
  890. # (row or column matrices), and z is a two-dimensional array. The
  891. # array z, is a function of x and y: z = f(x,y). Thus, the values in
  892. # the array x can be thought of a "row-labels", and the values of y
  893. # can be thought of as "column-lables" for the 2-dimensioal array z.
  894. #
  895. # At present plot3 can plot 3 distinct lists. Each list may have
  896. # different X, Y, and Z scales.
  897. #
  898. ##############################################################################
  899.  
  900. plot3 = function ( L31, L32, L33 )
  901. {
  902.   check_plot_object ();
  903.  
  904.   #
  905.   # 1st check list contents
  906.   #
  907.   
  908.   if (exist (L31)) { check_3d_list (L31); }
  909.   if (exist (L32)) { check_3d_list (L32); }
  910.   if (exist (L33)) { check_3d_list (L33); }
  911.   
  912.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  913.  
  914.   #
  915.   # Figure out the scale limits. 
  916.   # Needs improvement!
  917.   #
  918.   
  919.   xmin = xmax = ymin = ymax = zmin = zmax = 0;
  920.   if (exist (L31)) 
  921.   {
  922.     XYZ_scales (L31.x, L31.y, L31.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  923.     xmin = Xmin; xmax = Xmax;
  924.     ymin = Ymin; ymax = Ymax;
  925.     zmin = Zmin; zmax = Zmax;
  926.   }
  927.   if (exist (L32)) 
  928.   {
  929.     XYZ_scales (L32.x, L32.y, L32.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  930.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  931.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  932.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  933.   }   
  934.   if (exist (L33)) 
  935.   {
  936.     XYZ_scales (L33.x, L33.y, L33.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  937.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  938.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  939.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  940.   }
  941.   
  942.   _plgra ();
  943.   _plcol (1);
  944.   _pllsty (1);
  945.   _plfont (WIN.[P].font[p]);
  946.   _plwid (WIN.[P].width[p]);
  947.   
  948.   # basex = 2; basey = 2; height = 4;
  949.   xmin2d = -2.0; xmax2d = 2.0;
  950.   ymin2d = -3.0; ymax2d = 5.0;
  951.   
  952.   _plenv (xmin2d, xmax2d, ymin2d, ymax2d, 0, -2);
  953.   _plw3d (basex, basey, height, xmin, xmax, ymin, ymax, ...
  954.           zmin, zmax, WIN.[P].alt[p], WIN.[P].az[p]);
  955.   _plbox3 (WIN.[P].grid3x[p], WIN.[P].xlabel[p], 0, 0, ...
  956.            WIN.[P].grid3y[p], WIN.[P].ylabel[p], 0, 0, ...
  957.            WIN.[P].grid3z[p], WIN.[P].zlabel[p], 0, 0);
  958.   _plmtex ("t", 1.0, 0.5, 0.5, WIN.[P].title[p]);
  959.  
  960.   if (exist (L31))
  961.   {
  962.     if (find_char (WIN.[P].grid3x[p], "l"))
  963.     { x = log10 (real (L31.x)); else x = real (L31.x); }
  964.     if (find_char (WIN.[P].grid3y[p], "l"))
  965.     { y = log10 (real (L31.y)); else y = real (L31.y); }
  966.     if (find_char (WIN.[P].grid3z[p], "l"))
  967.     { z = log10 (real (L31.z)); else z = real (L31.z); }
  968.     
  969.     _plcol (2);
  970.     _plot3d (x, y, z, L31.x.n, L31.y.n, 3, 0);
  971.   }
  972.   if (exist (L32))
  973.   {
  974.     if (find_char (WIN.[P].grid3x[p], "l"))
  975.     { x = log10 (real (L32.x)); else x = real (L32.x); }
  976.     if (find_char (WIN.[P].grid3y[p], "l"))
  977.     { y = log10 (real (L32.y)); else y = real (L32.y); }
  978.     if (find_char (WIN.[P].grid3z[p], "l"))
  979.     { z = log10 (real (L32.z)); else z = real (L32.z); }
  980.     
  981.     _plcol (3);
  982.     _pllsty (2);
  983.     _plot3d (x, y, z, L32.x.n, L32.y.n, 3, 0);
  984.   }
  985.   if (exist (L33)) 
  986.   {
  987.     if (find_char (WIN.[P].grid3x[p], "l"))
  988.     { x = log10 (real (L33.x)); else x = real (L33.x); }
  989.     if (find_char (WIN.[P].grid3y[p], "l"))
  990.     { y = log10 (real (L33.y)); else y = real (L33.y); }
  991.     if (find_char (WIN.[P].grid3z[p], "l"))
  992.     { z = log10 (real (L33.z)); else z = real (L33.z); }
  993.     
  994.     _plcol (4);
  995.     _pllsty (3);
  996.     _plot3d (x, y, z, L33.x.n, L33.y.n, 3, 0);
  997.   }
  998.   
  999.   _plflush ();
  1000.   _pltext ();
  1001.   
  1002.   #
  1003.   # Increment the plot no. so that next time
  1004.   # we use the correct settings.
  1005.   #
  1006.   
  1007.   WIN.[P].subplot = WIN.[P].subplot + 1;
  1008.   
  1009.   return P;
  1010. };
  1011.  
  1012. ##############################################################################
  1013. #
  1014. # Plot a 3-Dimensional graph. The data is composed in a list, with
  1015. # elements `x', `y', and `z'. x and y are single-dimension arrays
  1016. # (row or column matrices), and z is a two-dimensional array. The
  1017. # array z, is a function of x and y: z = f(x,y). Thus, the values in
  1018. # the array x can be thought of a "row-labels", and the values of y
  1019. # can be thought of as "column-lables" for the 2-dimensioal array z.
  1020. #
  1021. # At present plmesh can plot 3 distinct lists. Each list may have
  1022. # different X, Y, and Z scales.
  1023. #
  1024. ##############################################################################
  1025.  
  1026. plmesh = function ( L31, L32, L33 )
  1027. {
  1028.   check_plot_object ();
  1029.   
  1030.   #
  1031.   # 1st check list contents
  1032.   #
  1033.   
  1034.   if (exist (L31)) { check_3d_list (L31); }
  1035.   if (exist (L32)) { check_3d_list (L32); }
  1036.   if (exist (L33)) { check_3d_list (L33); }
  1037.   
  1038.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1039.  
  1040.   #
  1041.   # Figure out the scale limits. 
  1042.   # Needs improvement!
  1043.   #
  1044.   
  1045.   xmin = xmax = ymin = ymax = zmin = zmax = 0;
  1046.   if (exist (L31)) 
  1047.   {
  1048.     XYZ_scales (L31.x, L31.y, L31.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  1049.     xmin = Xmin; xmax = Xmax;
  1050.     ymin = Ymin; ymax = Ymax;
  1051.     zmin = Zmin; zmax = Zmax;
  1052.   }
  1053.   if (exist (L32)) 
  1054.   {
  1055.     XYZ_scales (L32.x, L32.y, L32.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  1056.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  1057.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  1058.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  1059.   }   
  1060.   if (exist (L33)) 
  1061.   {
  1062.     XYZ_scales (L33.x, L33.y, L33.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  1063.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  1064.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  1065.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  1066.   }
  1067.  
  1068.   _plgra ();
  1069.   _plcol (1);
  1070.   _pllsty (1);
  1071.   _plfont (WIN.[P].font[p]);
  1072.   _plwid (WIN.[P].width[p]);
  1073.   
  1074.   # basex = 2; basey = 2; height = 4;
  1075.   xmin2d = -2.0; xmax2d = 2.0;
  1076.   ymin2d = -3.0; ymax2d = 5.0;
  1077.   
  1078.   _plenv (xmin2d, xmax2d, ymin2d, ymax2d, 0, -2);
  1079.   _plw3d (basex, basey, height, xmin, xmax, ymin, ymax, ...
  1080.           zmin, zmax, WIN.[P].alt[p], WIN.[P].az[p]);
  1081.   _plbox3 (WIN.[P].grid3x[p], WIN.[P].xlabel[p], 0, 0, ...
  1082.            WIN.[P].grid3y[p], WIN.[P].ylabel[p], 0, 0, ...
  1083.            WIN.[P].grid3z[p], WIN.[P].zlabel[p], 0, 0);
  1084.   _plmtex ("t", 1.0, 0.5, 0.5, WIN.[P].title[p]);
  1085.  
  1086.   if (exist (L31))
  1087.   {
  1088.     if (find_char (WIN.[P].grid3x[p], "l"))
  1089.     { x = log10 (real (L31.x)); else x = real (L31.x); }
  1090.     if (find_char (WIN.[P].grid3y[p], "l"))
  1091.     { y = log10 (real (L31.y)); else y = real (L31.y); }
  1092.     if (find_char (WIN.[P].grid3z[p], "l"))
  1093.     { z = log10 (real (L31.z)); else z = real (L31.z); }
  1094.     
  1095.     _plcol (2);
  1096.     _plmesh (x, y, z, L31.x.n, L31.y.n, 3);
  1097.   }
  1098.   if (exist (L32))
  1099.   {
  1100.     if (find_char (WIN.[P].grid3x[p], "l"))
  1101.     { x = log10 (real (L32.x)); else x = real (L32.x); }
  1102.     if (find_char (WIN.[P].grid3y[p], "l"))
  1103.     { y = log10 (real (L32.y)); else y = real (L32.y); }
  1104.     if (find_char (WIN.[P].grid3z[p], "l"))
  1105.     { z = log10 (real (L32.z)); else z = real (L32.z); }
  1106.     
  1107.     _plcol (3);
  1108.     _pllsty (2);
  1109.     _plmesh (x, y, z, L32.x.n, L32.y.n, 3);
  1110.   }
  1111.   if (exist (L33)) 
  1112.   {
  1113.     if (find_char (WIN.[P].grid3x[p], "l"))
  1114.     { x = log10 (real (L33.x)); else x = real (L33.x); }
  1115.     if (find_char (WIN.[P].grid3y[p], "l"))
  1116.     { y = log10 (real (L33.y)); else y = real (L33.y); }
  1117.     if (find_char (WIN.[P].grid3z[p], "l"))
  1118.     { z = log10 (real (L33.z)); else z = real (L33.z); }
  1119.     
  1120.     _plcol (4);
  1121.     _pllsty (3);
  1122.     _plmesh (x, y, z, L33.x.n, L33.y.n, 3);
  1123.   }
  1124.  
  1125.   _plflush ();
  1126.   _pltext ();
  1127.   
  1128.   #
  1129.   # Increment the plot no. so that next time
  1130.   # we use the correct settings.
  1131.   #
  1132.   
  1133.   WIN.[P].subplot = WIN.[P].subplot + 1;
  1134.   
  1135.   return P;
  1136. };
  1137.  
  1138. ##############################################################################
  1139. #
  1140. # Plot contours. The data is composed in a list, with
  1141. # elements `x', `y', and `z'. x and y are single-dimension arrays
  1142. # (row or column matrices), and z is a two-dimensional array. The
  1143. # array z, is a function of x and y: z = f(x,y). Thus, the values in
  1144. # the array x can be thought of a "row-labels", and the values of y
  1145. # can be thought of as "column-lables" for the 2-dimensioal array z.
  1146. #
  1147. ##############################################################################
  1148.  
  1149. plcont = function ( CL )
  1150. {
  1151.   check_plot_object ();
  1152.   
  1153.   #
  1154.   # 1st check list contents
  1155.   #
  1156.   
  1157.   if (exist (CL)) { check_3d_list (CL); }
  1158.  
  1159.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1160.   
  1161.   #
  1162.   # Figure out the scale limits. 
  1163.   # Needs improvement!
  1164.   #
  1165.  
  1166.   xmin = xmax = ymin = ymax = zmin = zmax = 0;
  1167.   if (exist (CL)) 
  1168.   {
  1169.     XYZ_scales (CL.x, CL.y, CL.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  1170.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  1171.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  1172.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  1173.   }
  1174.  
  1175.   _plgra ();
  1176.   _plcol (1);
  1177.   _pllsty (1);
  1178.   _plfont (WIN.[P].font[p]);
  1179.   _plwid (WIN.[P].width[p]);
  1180.  
  1181.   #
  1182.   # Set up the 1st viewport for drawing the plot.
  1183.   #
  1184.  
  1185.   if (!subplot_f) 
  1186.   {
  1187.     _pladv (0);        # Advance 1 subplot
  1188.   else
  1189.     subplot_f = 0;     # The user has set the subplot
  1190.   }
  1191.  
  1192.   _plvpas (0.15, 0.75, 0.15, 0.85, WIN.[P].aspect[p]);
  1193.   # _plvpor (0.15, 0.75, 0.15, 0.85);
  1194.   _plwind (xmin, xmax, ymin, ymax);
  1195.   _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  1196.  
  1197.   # Convert the data to log data if necessary.
  1198.   if (find_char (WIN.[P].gridx[p], "l"))
  1199.   { x = log10 (real (CL.x)); else x = real (CL.x); }
  1200.   if (find_char (WIN.[P].gridy[p], "l"))
  1201.   { y = log10 (real (CL.y)); else y = real (CL.y); }
  1202.   z = real (CL.z);
  1203.   
  1204.   if (exist (CL.clevel))
  1205.   {
  1206.     clevel = CL.clevel;
  1207.   else
  1208.     clevel = linspace(zmin, zmax, 10);
  1209.   }
  1210.  
  1211.   #
  1212.   # Draw the contours
  1213.   #
  1214.  
  1215.   l = 1;
  1216.   for (i in 1:clevel.n)
  1217.   {
  1218.     k = mod (i-1, 14) + 1;
  1219.     j = mod (i-1, 8) + 1;
  1220.     _pllsty(j);
  1221.     _plcol (1+k);
  1222.     if (_plcont (x, y, z, 1, CL.x.n, 1, CL.y.n, clevel[i]))
  1223.     {
  1224.       llevel[l] = clevel[i];
  1225.       l = l + 1;
  1226.     }
  1227.   }
  1228.  
  1229.   #
  1230.   # Reset color and draw the labels.
  1231.   #
  1232.  
  1233.   _plcol (1);
  1234.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  1235.  
  1236.   #
  1237.   # Draw the contour legend. Use a new viewport to the right
  1238.   # of the contour plot.
  1239.   #
  1240.  
  1241.   #_plvpas (0.75, 1.0, 0.15, 0.85, WIN.[P].aspect[p]);
  1242.   _plvpor (0.75, 1.0, 0.15, 0.85);
  1243.   _plwind (0, 1, 0, 1);
  1244.  
  1245.   v = 1 - 1/(2*llevel.n);
  1246.  
  1247.   for (i in 1:llevel.n)
  1248.   {
  1249.     xl = [0.1, 0.2, 0.3]';
  1250.     yl = [v, v, v]';
  1251.     v = v - 1/llevel.n;
  1252.  
  1253.     k = mod (i-1, 14) + 1;
  1254.     j = mod (i-1, 8) + 1;
  1255.  
  1256.     _plcol (1+k);
  1257.     _pllsty (j);
  1258.  
  1259.     _plline (3, xl, yl);
  1260.     sprintf (stmp, "%.2g", llevel[i]);
  1261.     plptex (stmp, xl[3]+.1, yl[3], , , 0);
  1262.   }
  1263.  
  1264.   # Flush  and go back to text mode.
  1265.   _plflush ();
  1266.   _pltext ();
  1267.   
  1268.   #
  1269.   # Increment the plot no. so that next time
  1270.   # we use the correct settings.
  1271.   #
  1272.   
  1273.   WIN.[P].subplot = WIN.[P].subplot + 1;
  1274.   
  1275.   return P;
  1276. };
  1277.  
  1278. ##############################################################################
  1279. #
  1280. # Plot 3-D lines, etc...
  1281. #
  1282.  
  1283. pl3d = function ( X, Y, Z, BR )
  1284. {
  1285.   local (X, Y, Z, BR)
  1286.   check_plot_object ();
  1287.   
  1288.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1289.  
  1290.   #
  1291.   # Some basic checks
  1292.   #
  1293.  
  1294.   if ((N = X.n) != Y.n) { error ("pl3d: X and Y must have same length"); }
  1295.   if (N != Z.n) { error ("pl3d: X and Z must have same length"); }
  1296.  
  1297.   if (!exist (BR)) { BR = N; }
  1298.   if (mod (N, BR) != 0) { error ("pl3d: X.n must be divisible by BR"); }
  1299.   iBR = int (N / BR);
  1300.   if (iBR == 1) { k = N; else k = BR; }
  1301.  
  1302.   #
  1303.   # Figure out the scale limits. 
  1304.   # Needs improvement!
  1305.   #
  1306.   
  1307.   xmin = xmax = ymin = ymax = zmin = zmax = 0;
  1308.   XYZ_scales (X, Y, Z, p, xmin, xmax, ymin, ymax, zmin, zmax);
  1309.  
  1310.   _plgra ();
  1311.   _plcol (1);
  1312.   _pllsty (1);
  1313.   _plfont (WIN.[P].font[p]);
  1314.   _plwid (WIN.[P].width[p]);
  1315.   
  1316.   if (!subplot_f) 
  1317.   {
  1318.     _pladv (0);        # Advance 1 subplot
  1319.     else
  1320.     subplot_f = 0;     # The user has set the subplot
  1321.   }
  1322.  
  1323.   if (find_char (WIN.[P].grid3x[p], "l"))
  1324.   { X = log10 (real (X)); else X = real (X); }
  1325.   if (find_char (WIN.[P].grid3y[p], "l"))
  1326.   { Y = log10 (real (Y)); else Y = real (Y); }
  1327.   if (find_char (WIN.[P].grid3z[p], "l"))
  1328.   { Z = log10 (real (Z)); else Z = real (Z); }
  1329.     
  1330.   # basex = 2; basey = 2; height = 4;
  1331.   xmin2d = -2.0; xmax2d = 2.0;
  1332.   ymin2d = -3.0; ymax2d = 5.0;
  1333.   
  1334.   _plvasp (WIN.[P].aspect[p]);
  1335.   _plwind (xmin2d, xmax2d, ymin2d, ymax2d);
  1336.   _plw3d (basex, basey, height, xmin, xmax, ymin, ymax, ...
  1337.           zmin, zmax, WIN.[P].alt[p], WIN.[P].az[p]);
  1338.   _plbox3 (WIN.[P].grid3x[p], WIN.[P].xlabel[p], 0, 0, ...
  1339.            WIN.[P].grid3y[p], WIN.[P].ylabel[p], 0, 0, ...
  1340.            WIN.[P].grid3z[p], WIN.[P].zlabel[p], 0, 0);
  1341.   _plmtex ("t", 1.0, 0.5, 0.5, WIN.[P].title[p]);
  1342.  
  1343.   _plcol (2);
  1344.   for (i in 1:iBR)
  1345.   {
  1346.     j = [(i-1)*k+1:i*k];
  1347.     _plline3 (k, X[j], Y[j], Z[j]);
  1348.   }
  1349.   _plflush ();
  1350.   _pltext ();
  1351.   
  1352.   #
  1353.   # Increment the plot no. so that next time
  1354.   # we use the correct settings.
  1355.   #
  1356.   
  1357.   WIN.[P].subplot = WIN.[P].subplot + 1;
  1358.   
  1359.   return P;
  1360. };
  1361.  
  1362. ##############################################################################
  1363. #
  1364. # error bar plot
  1365. #
  1366.  
  1367. plerry = function (x, y, y_low, y_high)
  1368. {
  1369.   check_plot_object ();
  1370.   
  1371.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1372.   WIN.[P].desc.[p] = 1j;    // use own legend
  1373.   
  1374.   if (x.nr != y.nr || x.nr != y_low.nr || x.nr != y_high.nr) 
  1375.   {
  1376.     error(" Size inconsistent in plerry.");
  1377.   }
  1378.   _plgra ();
  1379.   _plcol (1);
  1380.   _pllsty (1);
  1381.   _plfont (WIN.[P].font[p]);
  1382.   _plwid (WIN.[P].width[p]);
  1383.   xy_scales ( real([x,y,y_low,y_high]), p, xmin, xmax, ymin, ymax );
  1384.   
  1385.   if (!subplot_f) 
  1386.   {
  1387.     _pladv (0);        # Advance 1 subplot
  1388.   else
  1389.     subplot_f = 0;     # The user has set the subplot
  1390.   }
  1391.  
  1392.   if (WIN.[P].aspect[p] != 0)
  1393.   {
  1394.     _plvasp (WIN.[P].aspect[p]);
  1395.   else
  1396.     _plvsta ();
  1397.   }
  1398.  
  1399.   _plwind (xmin, xmax, ymin, ymax);
  1400.   _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  1401.  
  1402.   if (plot_matrix ( [x,y], p, 0, xmin, xmax, ymin, ymax, ymax-ymin ) < 0) 
  1403.   { 
  1404.     return -1; 
  1405.   }
  1406.  
  1407.   _plcol (3);
  1408.   _plerry(x.nr, x, y_low, y_high);
  1409.   _plcol (1);
  1410.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  1411.   _plflush ();
  1412.   _pltext ();
  1413.   
  1414.   #
  1415.   # Increment the plot no. so that next time
  1416.   # we use the correct settings.
  1417.   #
  1418.   
  1419.   WIN.[P].subplot = WIN.[P].subplot + 1;
  1420.   return P;      
  1421. };
  1422.  
  1423. ##############################################################################
  1424. #
  1425. # Plot a Histogram(s), from the columns of a matrix.
  1426. #
  1427. ##############################################################################
  1428.  
  1429. plhist = function ( M , nbin )
  1430. {
  1431.   check_plot_object ();
  1432.   
  1433.   if (!exist (nbin)) { nbin = 10; }
  1434.   
  1435.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1436.   np = M.nr;
  1437.   
  1438.   # Compute max/min values of data
  1439.   
  1440.   ymin = min (min (real (M)));
  1441.   ymax = max (max (real (M)));
  1442.   
  1443.   #
  1444.   # Check computed scale limits against user's
  1445.   #
  1446.   
  1447.   if (WIN.[P].ymin[p] != 1j) { ymin = WIN.[P].ymin[p]; }
  1448.   if (WIN.[P].ymax[p] != 1j) { ymax = WIN.[P].ymax[p]; }
  1449.  
  1450.   _plgra ();
  1451.   _plcol (1);
  1452.   _plfont (WIN.[P].font[p]);
  1453.   _plwid (WIN.[P].width[p]);
  1454.   
  1455.   for (i in 1:M.nc) 
  1456.   { 
  1457.     hscale[i] = hist_scales (M[;i], nbin);
  1458.   }  
  1459.  
  1460.   if (!subplot_f) {
  1461.     _pladv (0);        # Advance 1 subplot
  1462.   else
  1463.     subplot_f = 0;     # The user has set the subplot
  1464.   }
  1465.  
  1466.   if (WIN.[P].aspect[p] != 0)
  1467.   {
  1468.     _plvasp (WIN.[P].aspect[p]);
  1469.   else
  1470.     _plvsta ();
  1471.   }
  1472.  
  1473.   _plwind (ymin, ymax, 0, max (hscale));
  1474.   _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  1475.  
  1476.   v = max (hscale);
  1477.   xmax = ymax;
  1478.   for (i in 1:M.nc)
  1479.   {
  1480.     k = mod (i, 14) + 1;
  1481.     _plcol (WIN.[P].color[p;k]);
  1482.     _plhist (np, real(M[;i]), ymin, ymax, nbin, 1);
  1483.     
  1484.     if (!any (any (WIN.[P].desc.[p] == 1j)))
  1485.     {
  1486.       # Use the default if necessary
  1487.       if (WIN.[P].desc.[p][1] == "default") 
  1488.       {
  1489.     desc = "c"+num2str(i);
  1490.       else if (WIN.[P].desc.[p].n >= i) {
  1491.     desc = WIN.[P].desc.[p][i];
  1492.       else
  1493.     # Not sure what to do, user has messed up.
  1494.     desc = "";
  1495.       } }
  1496.  
  1497.       v = v - max(hscale)/11;
  1498.       xl = (ymax-ymin)*[10.5/12, 11/12, 11.5/12]' + ymin;
  1499.       yl = [v, v, v]';
  1500.  
  1501.       _plline (3, xl, yl);
  1502.       plptex(desc, xl[1]-(ymax-ymin)/25, yl[3], , , 1);
  1503.     }
  1504.   }
  1505.  
  1506.   _plcol (1);
  1507.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  1508.   _plflush ();
  1509.   _pltext ();
  1510.  
  1511.   #
  1512.   # Increment the plot no. so that next time
  1513.   # we use the correct settings.
  1514.   #
  1515.   
  1516.   WIN.[P].subplot = WIN.[P].subplot + 1;
  1517.  
  1518.   return 1;
  1519. };
  1520.  
  1521. ##############################################################################
  1522. #
  1523. # Various support functions for the WIN list
  1524. #
  1525. ##############################################################################
  1526.  
  1527. #
  1528. # Replot
  1529. #
  1530.  
  1531. replot = function ( )
  1532. {
  1533.   check_plot_object ();
  1534.   _replot ();
  1535. };
  1536.  
  1537. ##############################################################################
  1538. #
  1539. # Set the X-axis label
  1540. #
  1541.  
  1542. xlabel = function ( xstr )
  1543. {
  1544.   check_plot_object ();
  1545.   if (!exist (xstr)) { xstr = ""; }
  1546.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1547.   WIN.[P].xlabel[i] = xstr;
  1548. };
  1549.  
  1550. ##############################################################################
  1551. #
  1552. # Set the Y-axis label
  1553. #
  1554.  
  1555. ylabel = function ( xstr )
  1556. {
  1557.   check_plot_object ();
  1558.   if (!exist (xstr)) { xstr = ""; }
  1559.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1560.   WIN.[P].ylabel[i] = xstr;
  1561. };
  1562.  
  1563. ##############################################################################
  1564. #
  1565. # Set the Z-axis label
  1566. #
  1567.  
  1568. zlabel = function ( xstr )
  1569. {
  1570.   check_plot_object ();
  1571.   if (!exist (xstr)) { xstr = ""; }
  1572.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1573.   WIN.[P].zlabel[i] = xstr;
  1574. };
  1575.  
  1576. ##############################################################################
  1577. #
  1578. # Set the plot-title
  1579. #
  1580.  
  1581. pltitle = function ( xstr )
  1582. {
  1583.   check_plot_object ();
  1584.   if (!exist (xstr)) { xstr = ""; }
  1585.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1586.   WIN.[P].title[p] = xstr;
  1587. };
  1588.  
  1589. ##############################################################################
  1590. #
  1591. # Set the scale limits.
  1592. #
  1593.  
  1594. plimits = function ( xmin, xmax, ymin, ymax, zmin, zmax )
  1595. {
  1596.   check_plot_object ();
  1597.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1598.  
  1599.   if (exist (xmin)) 
  1600.   {
  1601.     WIN.[P].xmin[i] = xmin;
  1602.   else
  1603.     WIN.[P].xmin[i] = 1j;
  1604.   }
  1605.   if (exist (xmax)) 
  1606.   { 
  1607.     WIN.[P].xmax[i] = xmax;
  1608.   else
  1609.     WIN.[P].xmax[i] = 1j;
  1610.   }
  1611.  
  1612.   if (exist (ymin)) 
  1613.   {
  1614.     WIN.[P].ymin[i] = ymin;
  1615.   else
  1616.     WIN.[P].ymin[i] = 1j;
  1617.   }
  1618.   if (exist (ymax)) 
  1619.   {
  1620.     WIN.[P].ymax[i] = ymax;
  1621.   else
  1622.     WIN.[P].ymax[i] = 1j;
  1623.   }
  1624.  
  1625.   if (exist (zmin)) 
  1626.   {
  1627.     WIN.[P].zmin[i] = zmin;
  1628.   else
  1629.     WIN.[P].zmin[i] = 1j;
  1630.   }
  1631.   if (exist (zmax)) 
  1632.   {
  1633.     WIN.[P].zmax[i] = zmax;
  1634.   else
  1635.     WIN.[P].zmax[i] = 1j;
  1636.   }
  1637. };
  1638.  
  1639. ##############################################################################
  1640. #
  1641. # Set 2-D grid styles. A not-so-friendly interface.
  1642. #
  1643.  
  1644. plgrid = function ( sty_x, sty_y )
  1645. {
  1646.   check_plot_object ();
  1647.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1648.  
  1649.   if (exist (sty_x)) 
  1650.   { 
  1651.     if (class (sty_x) == "string")
  1652.     {
  1653.       WIN.[P].gridx[i] = sty_x;
  1654.     else
  1655.       error ("plgrid: requires string argument GRID_STY_X");
  1656.     }
  1657.   else
  1658.     WIN.[P].gridx[i] = grid_x_default;
  1659.   }
  1660.   if (exist (sty_y)) 
  1661.   { 
  1662.     if (class (sty_y) == "string")
  1663.     {
  1664.       WIN.[P].gridy[i] = sty_y;
  1665.     else
  1666.       error ("plgrid: requires string argument GRID_STY_Y");
  1667.     }
  1668.   else
  1669.     WIN.[P].gridy[i] = grid_y_default;
  1670.   }
  1671. };
  1672.  
  1673. ##############################################################################
  1674. #
  1675. # Set 3-D grid (axis) styles
  1676. #
  1677.  
  1678. plgrid3 = function ( sty_x, sty_y, sty_z )
  1679. {  
  1680.   check_plot_object ();
  1681.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1682.   if (exist (sty_x)) 
  1683.   { 
  1684.     if (class (sty_x) == "string")
  1685.     {
  1686.       WIN.[P].grid3x[i] = sty_x;
  1687.     else
  1688.       error ("plgrid3: requires string argument GRID_STY_X");
  1689.     }
  1690.   else
  1691.     WIN.[P].grid3x[i] = grid_3x_default;
  1692.   }
  1693.   if (exist (sty_y)) 
  1694.   { 
  1695.     if (class (sty_y) == "string")
  1696.     {
  1697.       WIN.[P].grid3y[i] = sty_y;
  1698.     else
  1699.       error ("plgrid3: requires string argument GRID_STY_Y");
  1700.     }
  1701.   else
  1702.     WIN.[P].grid3y[i] = grid_3y_default;
  1703.   }
  1704.   if (exist (sty_z)) 
  1705.   { 
  1706.     if (class (sty_z) == "string")
  1707.     {
  1708.       WIN.[P].grid3z[i] = sty_z;
  1709.     else
  1710.       error ("plgrid3: requires string argument GRID_STY_Z");
  1711.     }
  1712.   else
  1713.     WIN.[P].grid3z[i] = grid_3z_default;
  1714.   }
  1715. };
  1716.  
  1717. ##############################################################################
  1718. #
  1719. # A friendlier interface to changing 2-D grid/axis
  1720. # styles.
  1721. #
  1722.  
  1723. plaxis = function ( X_STR, Y_STR )
  1724. {
  1725.   check_plot_object ();
  1726.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1727.   
  1728.   if (exist (X_STR))
  1729.   {
  1730.     if (X_STR == "log") { WIN.[P].gridx[i] = "bcngstl"; }
  1731.   else
  1732.     WIN.[P].gridx[i] = grid_x_default;
  1733.   }
  1734.   
  1735.   if (exist (Y_STR))
  1736.   {
  1737.     if (Y_STR == "log") { WIN.[P].gridy[i] = "bcngstlv"; }
  1738.   else
  1739.     WIN.[P].gridy[i] = grid_y_default;
  1740.   }
  1741.   return P;
  1742. };
  1743.  
  1744. ##############################################################################
  1745. #
  1746. # Various internal support functions. Eventually these will be static.
  1747. #
  1748. ##############################################################################
  1749.  
  1750. #
  1751. # Find the X or Y scale limits .
  1752. # M can be a multi-column matrix, all columns
  1753. # will be used.
  1754. #
  1755.  
  1756. x_scales = function ( M, p, xmin, xmax )
  1757. {  
  1758.   #
  1759.   # 1st check for un-plottable data
  1760.   #
  1761.   
  1762.   if (any (any (isinf (M))))
  1763.   { error ("plot: cannot plot Infs"); }
  1764.   if (any (any (isnan (M))))
  1765.   { error ("plot: cannot plot NaNs"); }
  1766.   
  1767.   xmin = min (min (M));
  1768.   xmax = max (max (M));
  1769.   
  1770.   #
  1771.   # Check computed scale limits against user's
  1772.   #
  1773.   
  1774.   if (WIN.[P].xmin[p] != 1j) { xmin = WIN.[P].xmin[p]; }
  1775.   if (WIN.[P].xmax[p] != 1j) { xmax = WIN.[P].xmax[p]; }
  1776.  
  1777.   #
  1778.   # Check for potential errors
  1779.   #
  1780.   
  1781.   if (xmin == xmax) 
  1782.   { 
  1783.     # As good a guess as any
  1784.     xmin = xmin - 1;
  1785.     xmax = xmax + 1;
  1786.   }
  1787.   
  1788.   #
  1789.   # Finally, adjust if log-scales
  1790.   #
  1791.   
  1792.   if (find_char (WIN.[P].gridx[p], "l"))
  1793.   {
  1794.     if (xmin <= 0 || xmax <= 0) { error ("cannot plot log(x<=0)"); }
  1795.     xmin = log10 (xmin);
  1796.     xmax = log10 (xmax);
  1797.   }
  1798.   
  1799.   return 1;
  1800. };
  1801.  
  1802. y_scales = function ( M, p, xmin, xmax )
  1803. {
  1804.   
  1805.   #
  1806.   # 1st check for un-plottable data
  1807.   #
  1808.   
  1809.   if (any (any (isinf (M))))
  1810.   { error ("plot: cannot plot Infs"); }
  1811.   if (any (any (isnan (M))))
  1812.   { error ("plot: cannot plot NaNs"); }
  1813.   
  1814.   xmin = min (min (M));
  1815.   xmax = max (max (M));
  1816.   
  1817.   #
  1818.   # Check computed scale limits against user's
  1819.   #
  1820.   
  1821.   if (WIN.[P].ymin[p] != 1j) { xmin = WIN.[P].ymin[p]; }
  1822.   if (WIN.[P].ymax[p] != 1j) { xmax = WIN.[P].ymax[p]; }
  1823.  
  1824.   #
  1825.   # Check for potential errors
  1826.   #
  1827.   
  1828.   if (xmin == xmax) 
  1829.   { 
  1830.     # As good a guess as any
  1831.     xmin = xmin - 1;
  1832.     xmax = xmax + 1;
  1833.   }
  1834.   
  1835.   #
  1836.   # Finally, adjust if log-scales
  1837.   #
  1838.   
  1839.   if (find_char (WIN.[P].gridy[p], "l"))
  1840.   {
  1841.     if (xmin <= 0 || xmax <= 0) { error ("cannot plot log(y<=0)"); }
  1842.     xmin = log10 (xmin);
  1843.     xmax = log10 (xmax);
  1844.   }
  1845.   
  1846.   return 1;
  1847. };
  1848.  
  1849. z_scales = function ( M, p, xmin, xmax )
  1850. {
  1851.   
  1852.   #
  1853.   # 1st check for un-plottable data
  1854.   #
  1855.   
  1856.   if (any (any (isinf (M))))
  1857.   { error ("plot: cannot plot Infs"); }
  1858.   if (any (any (isnan (M))))
  1859.   { error ("plot: cannot plot NaNs"); }
  1860.   
  1861.   xmin = min (min (M));
  1862.   xmax = max (max (M));
  1863.   
  1864.   #
  1865.   # Check computed scale limits against user's
  1866.   #
  1867.   
  1868.   if (WIN.[P].zmin[p] != 1j) { xmin = WIN.[P].zmin[p]; }
  1869.   if (WIN.[P].zmax[p] != 1j) { xmax = WIN.[P].zmax[p]; }
  1870.  
  1871.   #
  1872.   # Check for potential errors
  1873.   #
  1874.   
  1875.   if (xmin == xmax) 
  1876.   { 
  1877.     # As good a guess as any
  1878.     xmin = xmin - 1;
  1879.     xmax = xmax + 1;
  1880.   }
  1881.   
  1882.   #
  1883.   # Finally, adjust if log-scales
  1884.   #
  1885.   
  1886.   if (find_char (WIN.[P].gridz[p], "l"))
  1887.   {
  1888.     if (xmin <= 0 || xmax <= 0) { error ("plot: cannot plot log(z<=0)"); }
  1889.     xmin = log10 (xmin);
  1890.     xmax = log10 (xmax);
  1891.   }
  1892.   
  1893.   return 1;
  1894. };
  1895.  
  1896. ##############################################################################
  1897. #
  1898. # Find the X and Y scales for a single matrix. (OLD)
  1899. #
  1900.  
  1901. xy_scales = function ( M, p, xmin, xmax, ymin, ymax )
  1902. {  
  1903.   #
  1904.   # 1st check for un-plottable data
  1905.   #
  1906.   
  1907.   if (any (any (isinf (M))))
  1908.   { error ("plot: cannot plot infs"); }
  1909.   if (any (any (isnan (M))))
  1910.   { error ("plot: cannot plot NaNs"); }
  1911.   
  1912.   if (M.nc == 1)
  1913.   {
  1914.     xmin = 1;
  1915.     xmax = M.nr;
  1916.     ymin = min (M);
  1917.     ymax = max (M);
  1918.   else
  1919.     xmin = min (M[;1]);
  1920.     xmax = max (M[;1]);
  1921.     ymin = min (min (M[;2:M.nc]));
  1922.     ymax = max (max (M[;2:M.nc]));
  1923.   }
  1924.   
  1925.   #
  1926.   # Check computed scale limits against user's
  1927.   #
  1928.   
  1929.   if (WIN.[P].xmin[p] != 1j) { xmin = WIN.[P].xmin[p]; }
  1930.   if (WIN.[P].xmax[p] != 1j) { xmax = WIN.[P].xmax[p]; }
  1931.   if (WIN.[P].ymin[p] != 1j) { ymin = WIN.[P].ymin[p]; }
  1932.   if (WIN.[P].ymax[p] != 1j) { ymax = WIN.[P].ymax[p]; }
  1933.   
  1934.   #
  1935.   # Check for potential errors
  1936.   #
  1937.   
  1938.   if (xmin == xmax) 
  1939.   { 
  1940.     xmin = xmin - 1;
  1941.     xmax = xmax + 1;
  1942.   }
  1943.   
  1944.   if (ymin == ymax)
  1945.   {
  1946.     ymin = ymin - 1;
  1947.     ymax = ymax + 1;
  1948.   }
  1949.   
  1950.   #
  1951.   # Finally, adjust if log-scales
  1952.   #
  1953.   
  1954.   if (find_char (WIN.[P].gridx[p], "l"))
  1955.   {
  1956.     if (xmin <= 0 || xmax <= 0) { error ("plot: cannot plot log <= 0"); }
  1957.     xmin = log10 (xmin);
  1958.     xmax = log10 (xmax);
  1959.   }
  1960.   if (find_char (WIN.[P].gridy[p], "l"))
  1961.   {
  1962.     if (ymin <= 0 || ymax <= 0) { error ("plot: cannot plot log <= 0"); }
  1963.     ymin = log10 (ymin);
  1964.     ymax = log10 (ymax);
  1965.   }
  1966.   
  1967.   return 1;
  1968. };
  1969.  
  1970. ##############################################################################
  1971. #
  1972. # Find the X, Y and Z scales for a single matrix.
  1973. #
  1974.  
  1975. XYZ_scales = function ( X, Y, Z, p, xmin, xmax, ymin, ymax, zmin, zmax )
  1976. {
  1977.   # X - scale
  1978.   if (any (any (isinf (X))))
  1979.   { error ("cannot plot infs"); }
  1980.   if (any (any (isnan (X))))
  1981.   { error ("cannot plot NaNs"); }
  1982.   
  1983.   xmin = min (real (X));
  1984.   xmax = max (real (X));
  1985.   
  1986.   # Y - scale
  1987.   if (any (any (isinf (Y))))
  1988.   { error ("cannot plot infs"); }
  1989.   if (any (any (isnan (Y))))
  1990.   { error ("cannot plot NaNs"); }
  1991.   
  1992.   ymin = min (real (Y));
  1993.   ymax = max (real (Y));
  1994.   
  1995.   # Z - scale
  1996.   if (any (any (isinf (Y))))
  1997.   { error ("cannot plot infs"); }
  1998.   if (any (any (isnan (Y))))
  1999.   { error ("cannot plot NaNs"); }
  2000.   
  2001.   zmin = min (min (real (Z)));
  2002.   zmax = max (max (real (Z)));
  2003.   
  2004.   #
  2005.   # Check computed scale limits against user's
  2006.   #
  2007.   
  2008.   if (WIN.[P].xmin[p] != 1j) { xmin = WIN.[P].xmin[p]; }
  2009.   if (WIN.[P].xmax[p] != 1j) { xmax = WIN.[P].xmax[p]; }
  2010.   if (WIN.[P].ymin[p] != 1j) { ymin = WIN.[P].ymin[p]; }
  2011.   if (WIN.[P].ymax[p] != 1j) { ymax = WIN.[P].ymax[p]; }
  2012.   if (WIN.[P].zmin[p] != 1j) { zmin = WIN.[P].zmin[p]; }
  2013.   if (WIN.[P].zmax[p] != 1j) { zmax = WIN.[P].zmax[p]; }
  2014.   
  2015.   #
  2016.   # Check for potential errors
  2017.   #
  2018.   
  2019.   if (xmin == xmax) 
  2020.   { 
  2021.     # As good a guess as any
  2022.     xmin = xmin - 1;
  2023.     xmax = xmax + 1;
  2024.   }
  2025.   
  2026.   if (ymin == ymax) 
  2027.   { 
  2028.     # As good a guess as any
  2029.     ymin = ymin - 1;
  2030.     ymax = ymax + 1;
  2031.   }
  2032.   
  2033.   if (zmin == zmax) 
  2034.   { 
  2035.     # As good a guess as any
  2036.     zmin = zmin - 1;
  2037.     zmax = zmax + 1;
  2038.   }
  2039.   
  2040.   #
  2041.   # Finally, adjust if log-scales
  2042.   #
  2043.   
  2044.   if (find_char (WIN.[P].grid3x[p], "l"))
  2045.   {
  2046.     if (xmin <= 0 || xmax <= 0) { error ("plot: cannot plot log(x<=0)"); }
  2047.     xmin = log10 (xmin);
  2048.     xmax = log10 (xmax);
  2049.   }
  2050.   
  2051.   if (find_char (WIN.[P].grid3y[p], "l"))
  2052.   {
  2053.     if (ymin <= 0 || ymax <= 0) { error ("plot: cannot plot log(y<=0)"); }
  2054.     ymin = log10 (ymin);
  2055.     ymax = log10 (ymax);
  2056.   }
  2057.   
  2058.   if (find_char (WIN.[P].grid3z[p], "l"))
  2059.   {
  2060.     if (zmin <= 0 || zmax <= 0) { error ("plot: cannot plot log(z<=0)"); }
  2061.     zmin = log10 (zmin);
  2062.     zmax = log10 (zmax);
  2063.   }
  2064.   
  2065.   return 1;
  2066. };
  2067.  
  2068. ##############################################################################
  2069. #
  2070. # Find the X and Y scales for a list of matrices
  2071. #
  2072.  
  2073. list_scales = function ( data, key, p, Xmin, Xmax, Ymin, Ymax )
  2074. {
  2075.   once = 1;
  2076.   
  2077.   for (i in members (data))
  2078.   {
  2079.     M = real (data.[i]);
  2080.     if (class (M) != "num") { continue; }
  2081.     
  2082.     if (abs (key) > M.nc)
  2083.     {
  2084.       error_1 ("plot: KEY argument > M.nc");
  2085.     }
  2086.     
  2087.     #
  2088.     # 1st check for un-plottable data
  2089.     #
  2090.     
  2091.     if (any (any (isinf (M))))
  2092.     { error ("plot: cannot plot infs"); }
  2093.     if (any (any (isnan (M))))
  2094.     { error ("plot: cannot plot NaNs"); }
  2095.     
  2096.     k = find ((1:M.nc) != abs (key));
  2097.     if (key > 0)
  2098.     {
  2099.       if (M.nc != 1)
  2100.       {
  2101.     x_scales ( real(M)[;key], p, xmin, xmax );
  2102.     y_scales ( real(M)[;k],   p, ymin, ymax );
  2103.       else
  2104.     x_scales ( (1:M.nr)', p, xmin, xmax );
  2105.     y_scales ( real(M),   p, ymin, ymax );
  2106.       }
  2107.     else if (key < 0) {
  2108.       x_scales ( real(M)[;k],        p, xmin, xmax );
  2109.       y_scales ( real(M)[;abs(key)], p, ymin, ymax );
  2110.     else
  2111.       x_scales ( (1:M.nr)', p, xmin, xmax );
  2112.       y_scales ( real(M),   p, ymin, ymax );
  2113.     } }
  2114.  
  2115.     if (once) 
  2116.     { 
  2117.       Xmin = xmin; Xmax = xmax; Ymin = ymin; Ymax = ymax; 
  2118.       once = 0; 
  2119.     }
  2120.     if (xmin < Xmin) { Xmin = xmin; }
  2121.     if (xmax > Xmax) { Xmax = xmax; }
  2122.     if (ymin < Ymin) { Ymin = ymin; }
  2123.     if (ymax > Ymax) { Ymax = ymax; }
  2124.   }
  2125.   
  2126.   return 1;
  2127. };
  2128.  
  2129. ##############################################################################
  2130. #
  2131. # Find the maximum number of elements in a bin for a single 
  2132. # column matrix.
  2133. #
  2134.  
  2135. hist_scales = function ( data, nbin )
  2136. {
  2137.   dmin = min (real (data));
  2138.   dmax = max (real (data));
  2139.   dbin = linspace (dmin, dmax, nbin+1);
  2140.   binval = zeros (nbin, 1);
  2141.   
  2142.   for (i in 1:nbin)
  2143.   {
  2144.     binval[i] = length (find (data >= dbin[i] && data < dbin[i+1]));
  2145.   }
  2146.   
  2147.   return max (binval);
  2148. };
  2149.  
  2150. ##############################################################################
  2151. #
  2152. # Plot the columns of a matrix (core function)
  2153. #
  2154. # Notes: This is the core function for plotting a matrix. If the
  2155. # matrix is a single column, then the matrix elements are plotted
  2156. # versus the row numbers. If it is a multi-column matrix, then
  2157. # columns 2:N are plotted versus column 1.
  2158. #
  2159. # p, K, k and l are indices for plot features.
  2160. # p: the current plot index (the plot #)
  2161. # K: usually 0. This index is used to start of the line style and
  2162. # color index (k = color index, l = line-style index). This is mostly
  2163. # used by plot_list, which may call plot_matrix repeatedly.
  2164. # k: the line color index. This value determines the line color used
  2165. # for each column of data. If K = 0, then k goes like 2:14, then
  2166. # flops back to 1:14.
  2167. # l: the line style inex. This value determines the line style used
  2168. # for each column of data - not the line-type (points, or lines). If
  2169. # K = 0, then l goes like 2:8, then flops back to 1:8.
  2170. #
  2171. ##############################################################################
  2172.  
  2173. plot_matrix = function ( M, key, p, K, xmin, xmax, ymin, ymax, v )
  2174. {
  2175.   np = M.nr;
  2176.   
  2177.   if (M.nc == 1)
  2178.   {
  2179.     x = 1:M.nr;
  2180.     y = real (M);
  2181.     k = mod (1+K, 14) + 1;
  2182.     l = mod (1+K, 8);
  2183.     
  2184.     if (find_char (WIN.[P].gridx[p], "l"))
  2185.     { x = log10 (x); }
  2186.     if (find_char (WIN.[P].gridy[p], "l"))
  2187.     { y = log10 (y); }
  2188.     
  2189.     _plcol (WIN.[P].color[p;k]);
  2190.     _pllsty (WIN.[P].lstyle[p;l]);
  2191.     
  2192.     if (get_style (WIN.[P].style.[p], k-1) == "line") 
  2193.     {
  2194.       _plline (M.nr, x, y);
  2195.     else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2196.       _plpoin (M.nr, x, y, WIN.[P].pstyle[p;l]);
  2197.     else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2198.       _plline (M.nr, x, y);
  2199.       _plpoin (M.nr, x, y, WIN.[P].pstyle[p;l]);
  2200.     else {
  2201.       _plline (M.nr, x, y);
  2202.     }}}}
  2203.  
  2204.     #
  2205.     # Now do the legend 
  2206.     #
  2207.     
  2208.     if (!any (any (WIN.[P].desc.[p] == 1j)))
  2209.     {
  2210.       # Use the default if necessary
  2211.       if (WIN.[P].desc.[p][1] == "default") 
  2212.       {
  2213.     desc = "c1";
  2214.       else if (WIN.[P].desc.[p].n >= k-1) {
  2215.     desc = WIN.[P].desc.[p][k-1];
  2216.       else
  2217.     # Not sure what to do, user has messed up.
  2218.     desc = "";
  2219.       } }
  2220.              
  2221.       v = v - (ymax-ymin)/11;
  2222.       xl = (xmax-xmin)*[10.5/12, 11/12, 11.5/12]' + xmin;
  2223.       yl = [v, v, v]' + ymin;
  2224.       
  2225.       if (get_style (WIN.[P].style.[p], k-1) == "line") 
  2226.       {
  2227.     _plline (3, xl, yl);
  2228.       else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2229.     _plpoin (3, xl, yl, WIN.[P].pstyle[p;l]);
  2230.       else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2231.     _plline (3, xl, yl);
  2232.     _plpoin (3, xl, yl, WIN.[P].pstyle[p;l]);
  2233.       } } }
  2234.  
  2235.       plptex(desc, xl[1]-(xmax-xmin)/25, yl[3], , , 1);
  2236.       
  2237.     }
  2238.  
  2239.   else
  2240.  
  2241.     #
  2242.     # Check for large column dimension
  2243.     #
  2244.     
  2245.     if (M.nc > 3*M.nr)
  2246.     {
  2247.       printf (" Plot %i columns and %i rows, are you sure (y/n) ? "...
  2248.                , M.nc, M.nr);
  2249.       ans = getline ("stdin");
  2250.       if (ans.[1] != "y") { return -1; }
  2251.     }
  2252.     
  2253.     ki = find ((1:M.nc) != abs (key));
  2254.     for (i in ki)
  2255.     {
  2256.       if (key > 0)
  2257.       {
  2258.     x = real (M[;key]);
  2259.     y = real (M[;i]);
  2260.       else if (key < 0) {
  2261.     x = real (M[;i]);
  2262.     y = real (M[;abs(key)]);
  2263.       else
  2264.     x = (1:M.nr)';
  2265.     y = real (M[;i]);
  2266.       } }
  2267.  
  2268.       # Check for log scales, adjust if necessary
  2269.       if (find_char (WIN.[P].gridx[p], "l"))
  2270.       { x = log10 (x); }
  2271.       if (find_char (WIN.[P].gridy[p], "l"))
  2272.       { y = log10 (y); }
  2273.       
  2274.       k = mod (i-1 + K, 14) + 1;
  2275.       l = mod (8 + i-2 + K, 8) + 1;
  2276.       
  2277.       _plcol (WIN.[P].color[p;k]);
  2278.       _pllsty (WIN.[P].lstyle[p;l]);
  2279.       
  2280.       if (get_style (WIN.[P].style.[p], k-1) == "line") 
  2281.       {
  2282.     _plline (np, x, y);
  2283.       else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2284.     _plpoin (np, x, y, WIN.[P].pstyle[p;l]);
  2285.       else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2286.     _plline (np, x, y);
  2287.     _plpoin (np, x, y, WIN.[P].pstyle[p;l]);
  2288.       else {
  2289.     _plline (np, x, y);
  2290.       }}}}
  2291.  
  2292.       #
  2293.       # Now do the legend 
  2294.       #
  2295.       
  2296.       if (!any (any (WIN.[P].desc.[p] == 1j)))
  2297.       {
  2298.     # Use the default if necessary
  2299.     if (WIN.[P].desc.[p][1] == "default") 
  2300.     {
  2301.       desc = "c" + num2str (i);
  2302.         else if (WIN.[P].desc.[p].n >= k-1) {
  2303.       desc = WIN.[P].desc.[p][k-1];
  2304.         else
  2305.       # Not sure what to do, user has messed up.
  2306.       desc = "";
  2307.         }}
  2308.              
  2309.     v = v - (ymax-ymin)/11;
  2310.     xl = (xmax-xmin)*[10.5/12, 11/12, 11.5/12]' + xmin;
  2311.     yl = [v, v, v]' + ymin;
  2312.     
  2313.     if (get_style (WIN.[P].style.[p], k-1) == "line") 
  2314.     {
  2315.       _plline (3, xl, yl);
  2316.         else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2317.       _plpoin (3, xl, yl, WIN.[P].pstyle[p;l]);
  2318.         else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2319.       _plline (3, xl, yl);
  2320.       _plpoin (3, xl, yl, WIN.[P].pstyle[p;l]);
  2321.         }}}
  2322.     
  2323.     plptex(desc, xl[1]-(xmax-xmin)/25, yl[3], , , 1);
  2324.     
  2325.       }
  2326.     }
  2327.   }
  2328.   
  2329.   return k-1;
  2330. };
  2331.  
  2332. ##############################################################################
  2333. #
  2334. # Plot all of the matrices in a list on the same plot
  2335. #
  2336.  
  2337. plot_list = function ( L, key, p, xmin, xmax, ymin, ymax )
  2338. {
  2339.   k = 0;
  2340.   v = ymax - ymin;
  2341.   
  2342.   #
  2343.   # Sort out the list members
  2344.   #
  2345.  
  2346.   sl = list_sort (L);
  2347.  
  2348.   # Plot the list members with numeric labels 1st.
  2349.   if (exist (sl.num))
  2350.   {
  2351.     for (i in sl.num)
  2352.     {
  2353.       M = L.[i];
  2354.       if (class (M) != "num") { continue; }
  2355.       if ((k = plot_matrix (M, key, p, k, xmin, xmax, ymin, ymax, v)) < 0) 
  2356.       { 
  2357.     return k; 
  2358.       }
  2359.     }
  2360.   }
  2361.  
  2362.   # Now plot the list members with string labels.
  2363.   if (exist (sl.char))
  2364.   {
  2365.     for (i in sl.char)
  2366.     {
  2367.       M = L.[i];
  2368.       if (class (M) != "num") { continue; }
  2369.       if ((k = plot_matrix (M, key, p, k, xmin, xmax, ymin, ymax, v)) < 0) 
  2370.       { 
  2371.     return k; 
  2372.       }
  2373.     }
  2374.   }
  2375.   return 1;
  2376. };
  2377.  
  2378. ##############################################################################
  2379. #
  2380. # Check the elements of LIST.
  2381. # LIST must contain elements `x', `y',
  2382. # and `z'
  2383. #
  2384.  
  2385. check_3d_list = function ( LIST )
  2386. {
  2387.   #
  2388.   # Check existence and types
  2389.   #
  2390.   
  2391.   if (class (LIST) != "list") {
  2392.     error ("plot3: argument must be a list");
  2393.   }
  2394.   if (!exist (LIST.x)) {
  2395.     error ("plot3: arg must contain `x' member");
  2396.   else if (class (LIST.x) != "num") {
  2397.     error ("plot3: x must be numeric");
  2398.   } }
  2399.   if (!exist (LIST.y)) {
  2400.     error ("plot3: arg must contain `y' member");
  2401.   else if (class (LIST.y) != "num") {
  2402.     error ("plot3: y must be numeric"); 
  2403.   } }
  2404.   if (!exist (LIST.z)) {
  2405.     error ("plot3: arg must contain `z' member");
  2406.   else if (class (LIST.z) != "num") {
  2407.     error ("plot3: z must be numeric");
  2408.   } }
  2409.  
  2410.   #
  2411.   # Check sizes
  2412.   #
  2413.   
  2414.   if (LIST.x.n != LIST.z.nr) 
  2415.   {
  2416.     error ("plot3: x.n != z.nr");
  2417.   }
  2418.   
  2419.   if (LIST.y.n != LIST.z.nc) 
  2420.   {
  2421.     error ("plot3: y.n != z.nc");
  2422.   }
  2423.   
  2424. };
  2425.  
  2426. ##############################################################################
  2427. #
  2428. # A special type of histogram plot.
  2429. #
  2430.  
  2431. plhistx = function ( M , nbin )
  2432. {
  2433.   check_plot_object ();
  2434.   
  2435.   if (!exist (nbin)) { nbin = 10; }
  2436.   
  2437.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  2438.   np = M.nr;
  2439.   
  2440.   # Compute max/min values of data
  2441.   
  2442.   ymin = min (min(real (M)));
  2443.   ymax = max (max(real (M)));
  2444.   
  2445.   #
  2446.   # Check computed scale limits against user's
  2447.   #
  2448.   
  2449.   if (WIN.[P].ymin[p] != 1j) { ymin = WIN.[P].ymin[p]; }
  2450.   if (WIN.[P].ymax[p] != 1j) { ymax = WIN.[P].ymax[p]; }
  2451.   
  2452.   _plgra ();
  2453.   _plcol (15);
  2454.   _plfont (WIN.[P].font[p]);
  2455.   _plwid (WIN.[P].width[p]);
  2456.   
  2457.   dbin = (linspace (ymin, ymax, nbin+1))';
  2458.   for (j in 1:M.nc) 
  2459.   {
  2460.     // counting
  2461.     for (i in 1:nbin) 
  2462.     {
  2463.       binval[i;j] = length (find (M[;j] >= dbin[i] && M[;j] < dbin[i+1]));
  2464.     }
  2465.   }
  2466.   
  2467.   if (!subplot_f) {
  2468.     _pladv (0);        # Advance 1 subplot
  2469.   else
  2470.     subplot_f = 0;     # The user has set the subplot
  2471.   }
  2472.  
  2473.   if (WIN.[P].aspect[p] != 0)
  2474.   {
  2475.     _plvasp (WIN.[P].aspect[p]);
  2476.   else
  2477.     _plvsta ();
  2478.   }
  2479.  
  2480.   xmin = 0;
  2481.   xmax =  max(max(binval));
  2482.   _plwind (ymin, ymax, xmin, xmax);
  2483.   _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  2484.  
  2485.   #
  2486.   # Now reorganize dbin and binval so we plot points in
  2487.   # the middle of the bins.
  2488.   #
  2489.  
  2490.   delbin = abs(dbin[2] - dbin[1]);
  2491.   dbin = (linspace (ymin+delbin/2, ymax-delbin/2, nbin))';
  2492.   dbin = [ymin ; dbin ; ymax];
  2493.   binval = [zeros(1,binval.nc); binval; zeros(1,binval.nc)];
  2494.  
  2495.   v = xmax;
  2496.   for (i in 1:M.nc)
  2497.   {
  2498.     k = mod (i, 14) + 1;
  2499.     l = mod (i,  8) + 1;
  2500.     _plcol (WIN.[P].color[p;k]);
  2501.     _pllsty (WIN.[P].lstyle[p;l]);
  2502.     
  2503.     if      (get_style (WIN.[P].style.[p], k-1) == "line") 
  2504.     {
  2505.       _plline (nbin+2, dbin, binval[;i]);
  2506.     else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2507.       _plpoin (nbin+2, dbin, binval[;i], WIN.[P].pstyle[p]+k);
  2508.     else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2509.       _plline (nbin+2, dbin, binval[;i]);
  2510.       _plpoin (nbin+2, dbin, binval[;i], WIN.[P].pstyle[p]+k);     
  2511.     } } }
  2512.  
  2513.     // write legend around upper-right corner.
  2514.     // it is better to have user to choose location for legend.
  2515.     
  2516.     if (!any (any (WIN.[P].desc.[p] == 1j)))
  2517.     {
  2518.       # Use the default if necessary
  2519.       if (WIN.[P].desc.[p][1] == "default") 
  2520.       {
  2521.     desc = "c"+num2str(i);
  2522.       else if (WIN.[P].desc.[p].n >= i) {
  2523.       desc = WIN.[P].desc.[p][i];
  2524.       else
  2525.     # Not sure what to do, user has messed up.
  2526.     desc = "";
  2527.       } }
  2528.  
  2529.       v = v - (xmax)/11;
  2530.       xt = (ymax-ymin)*[10.5/12, 11/12, 11.5/12]' + ymin;
  2531.       yt = [v, v, v]';
  2532.       
  2533.       if      (get_style (WIN.[P].style.[p], k-1) == "line") 
  2534.       {
  2535.     _plline (3, xt, yt);
  2536.       else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2537.     _plpoin (3, xt, yt, WIN.[P].pstyle[p]+k);
  2538.       else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2539.     _plline (3, xt, yt);
  2540.     _plpoin (3, xt, yt, WIN.[P].pstyle[p]+k);
  2541.       } } }
  2542.  
  2543.       plptex(desc, xt[1]-(ymax-ymin)/25, yt[3], , , 1);
  2544.     }
  2545.   }
  2546.   
  2547.   _plcol (15);
  2548.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  2549.   _plflush ();
  2550.   _pltext ();
  2551.   
  2552.   #
  2553.   # Increment the plot no. so that next time
  2554.   # we use the correct settings.
  2555.   #
  2556.   
  2557.   WIN.[P].subplot = WIN.[P].subplot + 1;
  2558.   
  2559.   return 1;
  2560. };
  2561.  
  2562.  
  2563. ##############################################################################
  2564. #
  2565. # Create a legend in the current plot window
  2566. #
  2567. # if pobj.desc.[p] = inf()        no legend
  2568. # if pobj.desc.[p] = "default"        default ("c1", "c2", ...)
  2569. # if pobj.desc.[p] = "string"        use "string" as description
  2570. #
  2571.  
  2572. #
  2573. # Set the current plot legend string
  2574. #
  2575.  
  2576. plegend = function ( LEGEND )
  2577. {
  2578.   check_plot_object ();
  2579.   
  2580.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  2581.   
  2582.   if (!exist (LEGEND)) 
  2583.   {
  2584.     WIN.[P].desc.[p] = 1j;
  2585.     return P;
  2586.   }
  2587.   
  2588.   if (class (LEGEND) == "string")
  2589.   {
  2590.     WIN.[P].desc.[p] = LEGEND;
  2591.   }
  2592.   
  2593.   return P;
  2594. };
  2595.  
  2596. set3d = function (bx, by, h)
  2597. {
  2598.   if (!exist (bx)) { basex = 2; else basex = bx; }
  2599.   if (!exist (by)) { basey = 2; else basey = by; }
  2600.   if (!exist (h)) { height = 4; else height = h; }
  2601. };
  2602.